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Abstract 

Quasi-stationary distributions, as discussed by Darroch & Seneta 
(1965), have been used in biology to describe the steady state be- 
haviour of population models which, while eventually certain to be- 
come extinct, nevertheless maintain an apparent stochastic equilib- 
rium for long periods. These distributions have some drawbacks: they 
need not exist, nor be unique, and their calculation can present prob- 
lems. In an earlier paper, we gave biologically plausible conditions 
under which the quasi-stationary distribution is unique, and can be 
closely approximated by distributions that are simple to compute. In 
this paper, we consider conditions under which the quasi-stationary 
distribution, if it exists, need not be unique, but an apparent stochas- 
tic equilibrium can nonetheless be identified and computed; we call 
such a distribution a quasi-equilibrium distribution. 
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1 Introduction 



A rather general population growth model can be formulated as a Markovian 
birth and death process X in continuous time, where X{t) represents the 
number of individuals at time t in a population in a prescribed area A, 
having transition rates 

,^ Qi^i+i = i/3{i/A), qi^i_i = i5{i/A), i > 1; 

Qij = otherwise, 

where and 6{x) are the per capita rates of birth and mortality at popu- 
lation density x = i/A. The logistic growth model of Verhulst (1838) was the 
first to describe mathematically the evolution of a population to a non-zero 
equilibrium, contrasting with the Malthusian law of exponential growth, and 
its stochastic version falls into the above framework, with = b constant 
in time, and with S{x) = d+cx, linearly increasing in x. However, if the set of 
states N is a communicating class and inf j,.>o S{x) > 0, the stochastic model 
does not have a non-zero equilibrium distribution even if /3(0) > 5(0) and 
lim^_j.oo(/3(a;) — S{x)) < 0, since N is then transient, and eventual absorption 
in is certain. 

Darroch and Seneta (1965), building on the work of Yaglom (1947) in the 
context of branching processes, introduced the concept of a quasi-stationary 
distribution, in an attempt to describe the long term behaviour of a tran- 
sient Markov chain prior to eventual absorption. However, for chains with 
countably infinite state space, Seneta and Vere- Jones (1966) showed that the 
quasi-stationary distribution need neither exist nor be unique. Furthermore, 
even when there is a unique quasi- stationary distribution, its calculation 
may pose substantial problems. This apparently makes the quasi-stationary 
distribution unsatisfactory for typical biological applications. However, in 
Barbour & PoUett (2010) [BP], we were able to give conditions, simply ex- 
pressed in terms of its properties, under which a continuous time Markov 
chain X has exactly one quasi-stationary distribution. Under the same con- 
ditions, the quasi-stationary distribution can be approximated to a specified 
accuracy by the equilibrium distribution of a positively recurrent 'returned 
process' X'^, which may often be much easier to compute. It was also shown, 
under slightly more stringent conditions, that the distribution of X{t) is close 
to its quasi- stationary distribution for long periods of time. 

The conditions given in [BP] are satisfied for many population models 
of the form f ll.ip . including that of Verhulst (1838). However, a related 
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model, in which the per capita death rate S{x) = d remains constant as x 
increases, and the birth rate dechnes exponentially, (3{x) = 6e~"^ for some 
a > (Ricker 1954), does not. Indeed, although this biologically plausible 
model also gives rise to apparently stable equilibrium behaviour for long 
periods of time, it follows from van Doorn (1991) that the process actually has 
infinitely many possible quasi- stationary distributions. To enable the long 
term behaviour of such models to be adequately described, we now introduce 
a new set of conditions, complementary to those in [BP], which can apply 
in cases, such as that above, in which the quasi-stationary distribution need 
not exist nor be unique. 

Denoting the state space of X by C U {0}, where is the cemetery state 
and C is irreducible, the returned process X^^ is also Markov. It evolves 
like X, except when it reaches the state 0. Whenever it does, instead of 
being absorbed in 0, it is instantaneously returned to C according to the 're- 
turn' probability distribution n; hence each X^ is a recurrent process. Under 
our conditions, the returned processes for a wide class of return distributions 
all have very similar equilibrium distributions, and the distribution of X{t), 
given any reasonable fixed initial state, is also similar to them for long peri- 
ods of time. Thus, for computational and practical purposes, the situation is 
much as before. The only difference is that the quasi- stationary distribution 
can no longer be taken as the representative of the class of 'good' equilib- 
rium distributions, since it need neither exist nor be unique. Instead, any 
member /i of the class of 'good' return distributions can be chosen, and the 
equilibrium distribution of then serves as a quasi-equilibrium approxima- 
tion to C{X{t)) in the appropriate range of t. 

The main results. Theorems 12.31 and 12. 4[ are proved in Section [2j In 
Section |31 as an illustration, we discuss the application of the theorems to 
birth and death processes. These processes have been widely studied, be- 
cause of their relatively simple structure, and allow our results to be easily 
interpreted. Our theorems are however equally applicable to processes with 
more complicated structure, and we illustrate their application to Markov 
population processes in several dimensions in Section |H 

2 The return approximation 

Assume that X is a stable, conservative and non-explosive pure jump Markov 
process on a countable state space, consisting of a single transient class C 
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together with a cemetery state 0. For any probability distribution /j, on C, 
define the modified process with state space C to have exactly the same 
behaviour as X while in C, but, on reaching 0, to be instantly returned to C 
according to the distribution /i. Thus, if Q denotes the infinitesimal matrix 
associated with X, and that belonging to X^, we have 

(2.1) q'^j = Qij + qiolij ioriJeC. 

In this section, under a rather simple set of conditions, we show that the 
stationary distribution vr^ of is little infiuenced by the choice of /i, for n in 
a large class Ai of distributions. We give a bound, uniform for all /x, G A4, 
on the total variation distance 

dTv(7r",7r^) := sup A} - n^iA}] = ^ J] l^'(^) ~ ^''(^)l 



Aec 



kec 



between ir'^ and n^, that is expressed in terms of hitting probabilities and 
mean hitting times for the process X. The bound is such that it can be 
expected to be small in circumstances in which the process X typically spends 
a long time in C in apparent equilibrium, before being absorbed in as a 
result of an 'exceptional' event. 
Define 

(2.2) TA := infO > 0: X{t) e A, X{s) i A for some s < t}, 

with the infimum over the empty set being taken to be oo, noting that > 
a.s. even when X{0) e A. Our basic conditions can then be expressed as 
follows. 

Condition B. There exists s e C such that, dehning 

Pk ^k[Xr^,^^y = s]; Tk := Efe[T{s,o}], 

we have 

(i) inf pfc = p > ; 

(ii) Tk < oo for all k eC. 

Here, and refer to the distribution of X conditional on X(0) — k. 

Condition B(i) can be expected to be satisfied in reasonable generality, 
and is the same as Condition A(i) in [BP]. Condition B(ii) substantially 
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relaxes Condition A(ii) in [BP], which stipulated that <T < oo, uniformly 
for all k E C. If X is typically to spend a long time in apparent equilibrium 
before being absorbed in 0, it will be necessary for 1 — ps, the probability 
that an excursion from s lands in 0, to be small. 
We first note that 

Ts = ^Tsk < oo, 

kec 

where 

poo 

Tsk / Fs[{T{s,o}>t}n{X{t)^k}]dt 
Jo 

is the expected amount of time spent in k before first returning to {s, 0}, 
starting in s. Hence, for any C > 0, we can pick C C such that 

(2.3) J]T,, < a^-Ps)Ts; 

we do so in such a way that s e C^, and that :— sup^g^ is as small 
as possible. We then define the process to be the same as X, except that 

any excursions outside take zero time to complete. This process now 
satisfies Condition A of [BP], so that the results of [BP] can be applied to it. 
Finally, we extend the results for X(^ to the process X. To accomplish this 
programme, we need some preparatory results. 

Lemma 2.1 Define fi(T) := J2kec f^(^)Tk- Then, under Condition B, 

(i) //(T) < Wt[^^ < i^iT)/p; 

(ii) Eferf^j < Tk + {l-pk)KT)/p, keC; 

(iii) IEst|*^| < oo if and only if ii{T) < oo , 

where is defined similarly to ta, hut with the process X^ in place of X, 
and E'* denotes expectation under the initial distribution /i. 

Proof. The proof is based on the equation 
(2-4) r^sy = + 

in which r^^jj is the time that elapses after r^^ until first reaches s, zero 
if X^{t{s,o}) — s. Taking expectations with respect to P'^, this yields 

= EM^)m+(i-Pfe)E%'f}, 

kec 
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from which it follows that 

/i(T) < E^rf^^ < ^(T) + (l-p)EVf;j 

and Part (i) is proved. Part (ii) follows by taking expectations in (12 ■4p with 
respect to P^, which also gives 

Efcrf^} > {l-pkWr^.y > (l-pfe)/i(r). 

Part (iii) follows from these considerations, taking k = s. ■ 

We now define 

Mm := {fieVMiCy. fi{T)<MT,}, 

for any M > 0. The next lemma bounds the equilibrium probability that 
^ C^, for any jj, e Mm- 

Lemma 2.2 Under Condition B, for any fi G Mm, we have 

vr'^iq) < {l-Ps){C + M/p} =: e{C,M). 

Proof. By a standard renewal argument, 

7r^(A)E,rf^^ = 5^ I r /[rf^^ > t] /[X(t) = k] dt\ 

k&A ^-^^ ^ 
keA 

It thus follows from ([23]), ([23D and Lemma EHi) that 

7r^{q)Ts < vr^(q)E.rf^} < (1 - p.)(CT. + /i(T)/p), 
and the lemma follows. ■ 



In what follows, we assume that M > 1, ensuring that the distribution Sg 
that puts probability 1 on the state s itself belongs to Mm- 
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We now return to the pure jump Markov process X^, which has the same 
jump chain as X, and the same jump rates for all k & Cq, but with 
Qk = oo for k ^ C(^. We also define its returned processes in the same 
way as for X, but with the new jump rates. We then define 

(2.5) TP := ^T,, > T4l-C(l-Ps)}, 

the mean time for to return to the set {0,s}, starting from s, the last 
inequality following from fl2.3p . 

Theorem 2.3 Suppose that Condition B holds, and that M > 1. Then, for 
any fi G Mm, 

ciTvK.Tr^O < 2(l-p,)f^ + C + -). 

\pTs P J 

Proof. We begin by considering the process Xj^ for any jj, ^ Mm-, noting 
that, for any k ^ C^, its equilibrium distribution vr^ satisfies 

(2.6) n^{k) = 7r^(A;)/7r^(C^). 

Now the process satisfies Condition A of [BP], and hence, from (2.13) 
of [BP], 

(2.7) rfTv(v^^^v^^'1 < 2(T+/p) ^ vrf (fc)gfco. 

Then, by a renewal argument, letting Nk{t) denote the number of visits of X''" 
to k in [0,t], we have 



y2n^'{k)qko = limt-iXo(t) 
keC 

(2.8) = \im{t~'Ns{t)}\im{No{t)/Nsit)} = Tr\l-ps). 

It now follows from i^M)-<^M that 

< 2(r+/p)T7i(l-p,)/7r^=(Cc). 
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Hence 

< 2(T+/p)T7i(l-p,) + 2£(C,M), 
this last from Lemma [221 as before, 5s G A^a/, because M > 1. I 



Remark. Of course, for the theorem to imply that tt^ and vr'^' are close, one 
needs (1 —p^) to be very small, which has already been noted as a necessary 

condition for long time stability. One also needs ^ + C + y not to be too 
large. The smaller C, is chosen, the larger is the value of , so that, in 
specific models, there is an optimum choice of limiting the accuracy of 
approximation that can be demonstrated by this method. 

We now turn our attention to the distribution of X{t) for fixed values 
of t, starting from any particular state in C^, and compare it to ■n^" . We 
begin by taking the initial state of X to be s, and remark later that this 
restriction makes little difference, provided that s is hit at least once. To 
state the theorem, we define 

:= Ws[X does not leave or hit before returning to s]. 

Since lim^_j.or^ = ps, the quantity 1 — can be made as close as desired to 
1 — Ps hj decreasing (, but at the cost of increasing at the same time. A 
crude bound for 1 — r.^ in terms of 1 — comes from observing that 

J^Tsk > (l-rc)/g(^), 

where q^'^^ := sup^gj^ Qk and 

:= {k ^ C^: q^j > for some j G C^}; 

from fl2.3p . this gives 
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Theorem 2.4 Suppose that Condition B holds, and let Bq := T^q^/p. If 
£(C, 1) < 1/2, then, for all t > 16T+/p, 

d^^{Cs{X{t)),7r'^) 
< (l-r^)(2t/T, + C + l/p) + DB^^J^+i2/er/'< =: vS)- 

Remark. Hence, informally, if (1 — rQ)B'^T^ / {pTs) ^ 1 and (1 — tq) ^ 1, 
the distribution Cs{X{t)) is close to tx^" for all times t such that 

BIT^Ip « t « T,/(l-rc); 

note that B(^ > 1, so that then t ^ /p also. 

Proof. The argument is based on coupling two copies Xf^ and Xf of the 

return process X^", with X^^^ in equilibrium and with X^^"* starting in s, 
by the method used in [BP], Theorem 2.5. The coupling is achieved by 
forcing X^^^ to follow the same sequence of states as X^^'^ after the first 
time that it hits s, and to have identical residence times in all states other 
than s; the careful matching of the exponentially distributed residence times 
of the two processes in s is all that is used to achieve the coupling. Now the 
argument leading to (2.18) of [BP] shows that X^^^ and X^^^ can be jointly 
defined in such a way that, if t > 16T^^/p, the event A^j that they have 
coupled by t is such that 

P[A^J < 4cGB^^+{2/er/''^t, 

for a universal constant cg, not depending on (. Now, because X^^-* is in 
equilibrium, 

P[X^^^ hits {0} U q before t] < tj^ ^ qki, 

k&Ci; ie{o}uc= 

and the double sum is bounded by (1 — rt^)/Ts''\ as in the argument leading 
to (EH). If X^^^ does not hit {0} U Q before t, and if A^t holds, then xf ^ 
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also avoids {0} U up to time t, in which case it is indistinguishable from 
an X-process, starting in s. It thus follows that 

(2.9) rf,v«%£.(X(t))) < (l-rc)(t/TP)+D5c^+(2/er/^<, 

with D := Acg- To complete the proof, it now merely remains to note that 
c?Tv(7r^% tt''') = 71^' (C^), and to use Lemma [521 note that (l-Ps) > (1 -''^c)' 
and that tP > TJ2, from ([23]), if e{(, 1) < 1/2. I 



Remark. Denoting by the event that X hits s before {0} U C^, the 
same argument can be used to show that (iTv('Cfc(X(t) | A^), £s(X(t))) is at 
most ?7((t) for any A; G C^, under the conditions of Theorem 12.41 Hence, con- 
ditional on the event that X hits s before reaching {0} UC^, the distribution 
of X(t) starting from any k E is also close to vr*^" for all times t such that 

B^J+/p « t « T,/(l-r^), 

provided that (1 — rc)B'^T^ j [pTg) ^ 1. Thus the return distribution ti^" is 
then indeed an appropriate long time approximation to the distribution of X 
in C, for times t -C Ts/(1 — r^), and vr'''' can be replaced by for any ix 
such that /i(T) < oo, with extra error at most that given by the bound in 
Theorem [isi with ^(T) /T, for M. 

The emphasis up to now has been on approximating C{X{t)) by tt*^". 
However, there are times when this approximation may also not be useful. 
Examples of this are processes in which a set can be found that has the 
properties that Tg^^ and T'^ are only moderately large and (1 — r^) is tiny, 
but for which (12.31) is not satisfied. Such is the case if there are states k ^ 
such that Tfc is extremely large; for instance, if the equilibrium around s is 
metastable, Tg itself may be enormously larger than T^'^K Here, nonetheless, 
the intermediate bound (12. 9|) shows that vr^"' acts as a good approximation 
for very long periods, even though vr''^ may be very different. 

In practice, computing vr^" may be complicated by having to cope with the 
detail of the return distribution from C^, which should not really be relevant 
here. The final approximation is therefore phrased instead in terms of the 
accelerated return process X^%for some C G C containing s but not 0, which 
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is returned to s at each time of leaving C. Here, the set C may reasonably be 
chosen to be finite, in which case computing the equilibrium distribution tt^, 
of the accelerated return process becomes relatively easy. We now define 
Tk,c' '■= ^ki'T'^ly], where f^" is defined as in f l2.2p . but with the process X^, in 

place of X; and we set T^, := sup f^^zc Tk,c' and tq' '■= ~ ^{s}]- 

Theorem 2.5 Suppose that Condition B (ii) holds, and let Be ■= T^'Qs- 
Then 



d,,{4^,C,{Xm < (l-fcO(t/r,,c') + ^5c'\M^ + (2/eWi6T+ 



T+ 

rv 



t 



for all t > 16TX, with D the same constant as in Theorem 2.4 



Proof. The argument runs exactly as in the proof of (12. 9p . but with the 
process X^J, instead of Since this process has no absorbing state 0, p 
can be replaced by 1 in the bound. I 



Theorem [23] is very much in line with the main message of the paper. The 
difference between Condition A of [BP] and Condition B of this paper largely 
concerns properties of the process starting from states that it rarely ever 
reaches, and such differences should not prevent effective approximation of 
the distribution of the process, at least for long periods of time. The essential 
difference between the situation in which Condition A is satisfied and that in 
which it is not is that, when it is not satisfied, the approximating distribution 
need not be a quasi-stationary distribution of the process, or even one of its 
return distributions, but instead a return distribution associated with the 
process restricted to a truncated state space. We consider an example of this 
in Section m 



3 Birth and death processes 

Let X be a birth and death process with birth rates bj > 0, 1 < j < oo, with 
6o = 0, and with strictly positive death rates dj, j > 0. Define ai = 1 and 
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then set 



Qm 



m -. 
^ OLldl 



l=r 



In order to use the theorems of the previous section, we need to find ex- 
pressions for the quantities p, ps, Tsk, Tg, and r,^ that appear there. These 
can be derived using hitting probabihties, which can be simply expressed us- 
ing the aj and the Sl^. First, for any j < m < I , we have 

(3.1) ¥m[X hits I before j] = S^+JS^+i. 

A first consequence is that 



(3.2) 
(3.3) 



1 -P^ 



P 



bs + ds 
Pi 



QS — l 



Sf 



1 



asips + ds)S{' 



diSf 



Next, if i ^ {0, s}, write Uki '■— ^k[T{i} < T{s,o}], ^ h "^u — 1- then we 
have 



r 

1 



(3.4) 



Uki = < 



if k < s < i or i < s < k; 
if s < i < k; 
if s < k < i; 



Sl/S\ if < A; < i < s; 

I Si^jSU, ifQ<i<k<s, 



from which it follows that, for such i, 

— TTZTT (W- + -Jrl if < i < s; 

ai{bi+di) \ 5?^^ Sl\ ' 



(3.5) 1 - P4t« < T{,,o}] 

Oii{hi+di)Slj^-^ 

Also, for i ^ {s, 0} and k s, we have 



if i > s. 



POO 

Tki := / ¥k[{r{s,o}>t}n{X{t)^i}]dt 
Jo 



(3.6) 



Uki 



bi + di 1 -Pi[T{i} < T{s,o}]' 
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with Tfcs = Tfco = 0, and then T^g = l/{bs + dg) and 

T^s — 1 , i 



(3.7) T,, 
These in turn give 

(3.8) Tk = Y,Tki = 



bs+ds 
bs+ds 



if < 2 < s; 
if i > s. 



i>i [ J2i>s+i^iSi+i if /OS, 

and 

(3.9) T, = (6,T,+i + 4T,_i + l)/(6, + 4) = \ ^ . y^a^SY^ySl 

Now, choosing any value of s > 0, the formulae (13. 7p . (13. 2 p and (13. 9 p can be 
used for any ( to determine a suitable set := {1,2,..., a^} so that (12. 3p is 
satisfied, and (13. 8p enables both and fi(T) to be determined. Furthermore, 
it follows from (13. 5 p with for s and with s for i that 

Thus, and from (13. 3p . all the ingredients for the bounds in Theorems 12.31 
and l2.4l are available, recalling also, for the calculation of S^, that Qs = bs+dg. 

For example, take the birth and death process given in (II. ip with A large, 
S{x) = d constant, and with /?(■) given by the Ricker choice (3{x) = be~°'^; 
thus bj = jbe~°'^^^ and dj = jd. If b > d, the deterministic equilibrium, in 
which the birth and death rates are equal, is given by x = ^ \og{b/d) =: c > 0, 
suggesting the choice of s := s{A) := [Ac\. This gives 



a 



lg-Qi/2Ag-aO-s)V2A. 



(3.11) aj/{sas} X j' 



note that is exponentially large in A. Thus, immediately, > 1/d, and 
> l/d if a > 2{s{A) + 1). Hence I - Ps = 0{{b/dy'/'^} from (Q, and 
l-r^ = 0{{b/d)-'/^} also if > 2{s{A) + 1), from fl3^ . Furthermore, S( 
is uniformly bounded in A, so that p is uniformly bounded below, by (13. 3p . 
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To choose the set := {1,2,..., a^}, note that, from fl3.4p -( 1X7|) and fl3.1ip . 

a sips + ds) 



for all i, and thus 



i>l 

Hence, from (12. 3p . the choice [a^ = 2(s + 1)] corresponds to a value of C < 1- 
For the corresponding value of T^, it is necessary to bound the expressions 
for Tfc, which is in detail tedious; however, it is not difficult to deduce that 
Tfc = 0(logA V log/c), so that = O(logA). From Theorem 12. 4[ it now 
follows that £s(X(t)) is close to vr''^ for all times t such that 

A\logAf < t < A-^/2gA/2a_ 

Furthermore, fi{T) = 0(logy4) for all return distributions concentrated on 
sets of the form {1, 2, . . . , A"^} for any fixed exponent m, and Theorem 12.31 
thus shows that the corresponding equilibrium distributions vr^ are all ex- 
ponentially close to TT*^' as A — )■ oo — indeed, /z would have to have ex- 
traordinarily long tails for anything else to be the case. Hence the fact that 
this process has infinitely many quasi-stationary distributions should not be 
interpreted as showing any kind of practical instability, at least for large A: 
there is an effective long time stable distribution, and it is extremely close 

to 7?"^=. 

Rather similar analyses could be undertaken for a variety of other well- 
known models. An analogue of the Beverton & Holt (1957) model would 
have f3{x) = 6/(1 + x/m) for 6 > and m > 0, that of Hassell (1975) would 
have P{x) = 6/(1 + x/mY, and that of Maynard-Smith & Slatkin (1973) 
would have (3{x) = 6/(1 + (x/mY). The qualitative conclusions would be 
entirely similar. 



4 Markov population processes 

In this section, we consider Markov population processes := {X^it), t > 
0), > 1, taking values in Z'^, for some d > 1. In many applications, 
the components represent the numbers of individuals of a particular type or 
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species, with a total of d types possible. The process evolves as a Markov 
process with state-dependent transitions 

(4.1) X -> X + J at rate Naj{N-^X), X e Z^f , J e J, 

where JT" C Z'^ is a fixed finite set, and we define J* := maxjgj- 1 J|. Density 
dependence is refiected in the fact that the arguments of the functions aj are 
counts normalised by the 'typical size' N. The functions aj: Z4. R+ are 
assumed to be twice continuously differentiable on M^j_, and to be such as to 
ensure that Xjv is locally irreducible; that is, the number of steps required to 
get from any state X 7^ to any of its lattice neighbours X + e^^\ 1 < j < (i, 
with positive probability, is uniformly bounded. 

Such processes satisfy a law of large numbers (Kurtz, 1970), expressed in 
terms of the system of deterministic equations 

(4.2) f = 5^J«,(0 =: F{0, ^eR'; 

here, ^{t) approximates XAr(t) := N^^X]^{t), and the quantity F represents 
the infinitesimal average drift of the components of the random process. We 
now suppose that F(c) = for some c G M*^ with Cj > 0, 1 < j < (i, and that 
all the eigenvalues of the matrix of derivatives DF{c) =: A have negative 
real parts. In this case, c is a locally stable equilibrium of the deterministic 
system (14. 2p . and, if is started with X~^XAr(0) close to c, the law of large 
numbers implies that XAr(t) remains close to c, in the sense that 

(4.3) sup \xiy(t) — c| -^d 0, 

0<t<T 

for any finite T > 0. The central limit theorem in Kurtz (1971) also shows 
that 

(4.4) N^/'^{xn{-) - c) ^ X in D[0,T], 

for any T > 0, where x is a Gaussian process whose stationary distribution 
has zero mean and covariance matrix S satisfying 

(4.5) AS + SA^ + ct2(c) = 0, 

where a'^{x) := J^jej'^'^^^J^^)- Here, we complement this approximation, 
by using Theorem 12.51 to show that the distribution of X^it) is close in total 
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variation, for time periods that become extremely long as increases, to the 
equilibrium distribution vr^ of a truncated process Xn, which is returned to 
a specified state s := sn near Nc whenever it leaves a neighbourhood C'{N) 
of Nc. Of course, this distribution, appropriately centred and normalized, 
converges to MVNd{0, S) as — > oo. 

In order to prove such a result, we need to define the neighbourhood C'{N), 
and to show that the quantities (1 — rc'{N)), '^/Tc'{n) and T^i^^j^-^ appearing 
in Theorem 12.51 can be suitably bounded. The inequality 

(4.6) l/fc'iN) < Qsr, = nJ^^^AN-'sn) 

is immediate. For the remaining bounds, we use Lyapounov-Foster-Tweedie 
methods (Meyn & Tweedie, 1993). We write y := N^^'^{x — c) and := 
y'^Vy, where the positive definite symmetric matrix V is to be chosen later, 
and we first consider the process 1/Ar(-) := N^^'^{xn{-) — c) stopped when 
kyjvWl — CqA^^''^, for Co also to be chosen later. Then, for y such that \ry\ < 
cqN^^"^, the generator A of the Markov process acting on a real function g{y) 
takes the form 



(4.7) {Ag){y) = Najic + N-'/'y){giy + N-'/'j) - g{y)}. 
Using Taylor's expansion on ^f, for |r| < CqA^^^^, we have 

{Ag){y) - 5^iVa^(c + iV-V2^){iV-V2j^Z)^7(l/) + \N-^D'g{y)[J^% 

(4.8) < N-^l%{y-g), 
where 

r/3(y;^7) := ^afjl sup \\D^g{y + u)\l 

jej |h|<JV-i/2j* 

and a}° := sup^.(2,_c)Ty(^_c)<c2 a;j(x). Similarly, expanding aj, we obtain 



(4.9) 



iVV2 ^ aj{c + N-~'"y)J^Dg{y) - y^A^Dg{y) 
JeJ 



< N-'/%iy;g), 
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where we have used the facts that Xljej ■Joiji'^) = -^(c) = and that 
J2jp7 = where 



r^,{y;g) := J, 



and af := sup^^(^_^)Ty(^_,)<c2 and then 



(4.10) 
with 



5^{aj(c + N-'/^y) - ajic)}D^g{y)[J, J] 
JeJ 



< N-'/%{y;g), 



V2{y;9) ■■= \\D'giy)\\J',J2'^j\y\ 



JeJ 

and a*/ := sup2..(^_c)Ty(^_c)<c2 ||-Dq;j(x)||. Thus, if one ignores the error 
terms, the generator acts on g as that of a multivariate Ornstein-Uhlenbeck 
process, 

(4.11) {Ag){y) ^ y'^A^Dgiy) + ^tr {a {c)D^g{y) a (c)}, 

with drift matrix A and infinitesimal covariance matrix cr^(c). 

We now consider the generator acting on functions g of the form g{y) = 
Fe{r), where F^ij) := f(t)dt, and the function / is non-negative. This 
gives 

Dg{y) = f{r)r;'Vy; D'g{y) = r;'f{r)V+{f'iry)-r;'f{ry)}rfVyiVyf 

Thus the first term in the approximation (14. lip to A yields 

y^A^Dgiy) = firy)r;V A'^Vy = lf{ry)ryV{A^V + VA}y. 

In order to choose functions g such that g{yNit)) is a super-martingale, we 
would like the right hand side to be negative, which will be the case if V 
is chosen in such a way that the symmetric matrix [A'^V + VA) is negative 
definite. One way of doing so here is to take V := where S is as in (14. 5p . 
in which case 

A^V + VA = S-^{SA^ + AS}S-^ = -S-V2(c)S-^ 
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is immediately negative definite. The remaining term in f l4.1ip then gives 
\ii{a{c)D'g{y)a{c)} 



= ir^-V(r,)tr{a(c)l-a(c)} 

+ H/'(^.) V(r,)K-2tr{a(c)Vi/(Vi/)^a(c)} 
= i^;V(r,)tr {a'} + - r^' firy)}Riy) , 

where 5"^ := S~^/^(T^(c)S~^''^ is positive definite, and 



R{y) 



y'^T. ^y 



is bounded between its smallest and largest eigenvalues 7 and T. 

We begin by taking /(r) := r'^^e^^ , for m and /3 to be chosen suitably. 
Then (14. lip gives the main part of {Ag){y) as 



2^ ' y 



{-r^V^-'^'^'E-'^'y + r-hr {a^ + {2/3 - r;'{m + l)}R{y) } 

< -ie^'-^r;™+i{-(7 - 2/3) - ry\m + 1 - tr {a^})} 
=: -G(r,), 

say. We now choose /3 and m in such a way that 2/3 < 7 and m + 1 > tr {5"^}. 

For the remainders, we note first, for rj^, that there exist constants c_, c+ 
and i^' such that 

sup \\D^g{y + u)\\ < K\\D^g{y)\\ for c.iV"^/^ < ||y|| < c+A^^/^ 

|M|<Ar-i/2j, 

and that HD^^fd/)!! < Csr^ exp(/3r^){l + r^} for some C3. Thus, for 
all y such that N~^/'^cs < ry < c'^N^/"^, for suitable 03,03, where we also 
choose C3 < Co, it follows that rj^^y; g) < G{ry)/6. Similar considerations for 
rii{y] g) and rj2{y;g) show that, possibly increasing C3 and decreasing C3, the 
inequality 

Vi{y; 9) + V2{y; 9) + risiy; g) < lG{ry) 
holds for all y G fi(A^~^/2c3, where 

B{p,R) := {y: p < r, < i?}. 

Hence, for such y, we always have 

{Ag){y) < -lG{ry)<0. 



Thus the quantity F^iy^it A fpRj) is a non- negative super-martingale, for 

any A^~^/^C3 < p < R < N^^'^c'^ and any < £ < p, where 

V« := inf{y^(t)^5(p,i?)}. 

Defining p{p,R]r) := P[|/Ar(v/j) G 5(0, p) | ry^(o) = r], it thus follows 
easily from the optional stopping theorem that 



m+l 



4/9 /?' 

(4.12) l-p(P,^;0 < F,{r)/F,{R) < ^-^^^ e 



for p, i? such that N-^/^c^ <p<R< N^^'^c'^ and 2(3 R{R - p) > 1, with the 
last condition ensuring that a simple lower bound for Fp{R) is valid. So take 

(4.13) C'(iV) := {X:N-'^\X - Nc) e B{0,N'/^c'^)}, 

and let stv € Z^^ be the closest lattice point to Nc. Any path of ?/7v starting in 
5(0, A^^^/^Ca) has positive probability of hitting N~^/'^{sn — Nc) by taking 
the most direct path from ?/Ar(0) to N~^/'^(sn — Nc), and this probability 
is uniformly bounded away from 0, by the local irreducibility assumption 
on Xn, and because the number of possible values of yN{-) in -8(0, N'^^'^c^) 
is uniformly bounded as N varies. Hence 

(4.14) P[r|,^| < fo,2iv-V2e3 I VNiO) G 5(0, N'^/'c^)] > 6 

for some 6 > 0. If the complementary event occurs, then 5(0, N'^^'^c^) is hit 
again by j/at before it leaves 5(0, A^^/^Cg) with probability at least 

1 - KN"" e-(^3)'/3^, 
for some K, in view of fl4.12p . It thus follows that 

(4.15) l-rc'iN) < r ^iTAT'- e-(^3m = 0(6"'^'^), 

for any < /3' < (4)^/3. 

In order to control the mean time to hitting s for the process X^v, we 
take /(r) := r"*" + Or, for m large enough and 6 small enough positive. 
Then {Ag){y) once again has two principal negative contributions, the first, 
bounded above by — ^T^'r^, coming from the drift term, and the second. 
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bounded above by —^mvy , from the variance term. The former domi- 
nates all positive terms for Vy > vq, for some fixed tq, and the second then 
dominates for the smaller values of Vy, if m is chosen large enough; the quan- 
tities rii{y;g), 1 < / < 3, are treated as before, and the upper and lower 
bounds for ry can be left unchanged. Hence, in the same range of y, we 
always have 

{Ag)iy) < -S'<0, 
for some 5' > 0. Applying the optional stopping theorem then yields 

m%,n\ry^io)=r} < F,{r) < ^ + 

if m > 1, uniformly in p, R such that N'^^'^cs < p < R< N^/'^c'^. Take the 
extreme values for p and R. Then since, for this i?, the process is returned 
directly to sat if yN{jp,R) ^ -8(0, i?), and since the mean time to either hitting 
N-^/'^{sn - Nc) or leaving B{{],2N-^/^c^), starting within 5(0, A^-^/^ca), is 
uniformly bounded by some Ci < cxd, a regenerative argument much as above 
shows that 

< co{Ar(-^)/2 + Ar} + ci-t-(l-5)f+(^), 
for 5 as in fl4.14p . and hence that, uniformly in A^, 
(4.16) f+(^) < C{Ar(™-i)/2 + Ar}, 

for some C < oo. 

Collecting the above bounds, we have enough to prove the following the- 
orem. 

Theorem 4.1 Suppose that X is a Markov population process with tran- 
sition rates Naj as given in (14.11) . and that the aj are such as to ensure 
that Xn is locally irreducible. Suppose also that F[c) = for some c G M"^ 
with Cj > 0, I < j < d, and that all the eigenvalues of the matrix of deriva- 
tives DF{c) have negative real parts. Then there exist a, f3i, (32 and c'^ > 0, 
and Ci, C2 and C3 < 00, depending only on the parameters of the process 
and not on N, such that, for all t, 

d^,{n'^,^^^,Cs{Xm < C,tNe-^^'' + C2t-'/'N'^-'"/' + Cse-P^'''-\ 

where s = is the nearest lattice point to Nc, and C'{N) is as defined 
%n 
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Proof. All that is needed is to apply the estimate given in Theorem 12.51 An 
upper bound on (1 — rc'{N)) is given in fl4.15l) : a bound on is given 

in fl4.16p : and 1/Tc'(n) is bounded in fl4.6p . For the exponent f3i, any (3' as 
for fl4.15p can be taken; a = max{(m — l)/2, 1} as in fl4.16p : and (32 can be 
taken to be (1 - log2)/{32C}, for C as in KW} . ■ 



In view of Theorem 14.11 the equilibrium distribution vr^ is a very good 
approximation in total variation to £s(Xjv(t)), provided that t is bounded be- 
low by a suitable power of N and above by a quantity growing exponentially 
with N. 

The lower bound given here for the time at which the quasi-equilibrium 
approximation becomes accurate is very pessimistic. The main reason is 
that the general coupling strategy used to prove Theorems 12.41 and 12.51 can 
be very inefficient in specific instances, and is so here. Better results could be 
expected by using the methods to be found in Roberts & Rosenthal (1996). 
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